library(tidyverse)
library(knitr)pal_control = c("#2E86B4", "#F2B827")
pal_outcome = wesanderson::wes_palette("Zissou1", 3, "continuous")
pal_condition = wesanderson::wes_palette("Zissou1", 6, "continuous")
pal_condition8 = wesanderson::wes_palette("Zissou1", 16, "continuous")
pal_condition8 = c(pal_condition8[1], pal_condition8[3],
pal_condition8[14], pal_condition8[5],
pal_condition8[15], pal_condition8[16],
pal_condition8[9], pal_condition8[12])source("load_data.R")
ind_diffs = ind_diffs %>%
mutate(bmi_z = as.numeric(scale(bmi)),
fat_z = as.numeric(scale(fat)))# roi distributions
betas %>%
filter(session == "all") %>%
select(subjectID, con, condition, control, roi, process, meanPE) %>%
unique() %>%
filter(!is.na(roi)) %>%
ggplot(aes(meanPE, fill = roi)) +
geom_density(color = NA) +
facet_wrap(~roi, ncol = 4) +
theme_minimal(base_size = 14) +
theme(legend.position = "none")# dot product distribution
dots %>%
filter(session == "all") %>%
select(subjectID, con, condition, control, process, test, dotProduct) %>%
unique() %>%
ggplot(aes(dotProduct, fill = test)) +
geom_density(color = NA) +
facet_wrap(~process, ncol = 3, scales = "free") +
theme_minimal(base_size = 14) +
theme(legend.position = c(.85, .2))dots %>%
filter(session == "all") %>%
select(subjectID, con, condition, control, process, test, dotProduct, mask) %>%
unique() %>%
ggplot(aes(dotProduct, fill = mask, alpha = .5)) +
geom_density(color = NA) +
facet_wrap(~process, ncol = 3, scales = "free") +
theme_minimal(base_size = 14) +
theme(legend.position = c(.85, .2))# standardize and winsorize
betas_std = betas %>%
group_by(roi, session) %>%
mutate(meanPE_std = scale(meanPE, center = TRUE, scale = TRUE),
outlier = ifelse(meanPE_std > 3, "yes",
ifelse(meanPE_std < -3, "yes", "no")),
meanPE_std = ifelse(meanPE_std > 3, 3,
ifelse(meanPE_std < -3, -3, meanPE_std))) %>%
ungroup()
dots_std = dots %>%
group_by(map, test, mask, session) %>%
mutate(dotProduct_std = scale(dotProduct, center = TRUE, scale = TRUE),
outlier = ifelse(dotProduct_std > 3, "yes",
ifelse(dotProduct_std < -3, "yes", "no")),
dotProduct_std = ifelse(dotProduct_std > 3, 3,
ifelse(dotProduct_std < -3, -3, dotProduct_std))) %>%
ungroup()# summarize
betas_std %>%
filter(session == "all") %>%
group_by(outlier) %>%
summarize(n = n()) %>%
spread(outlier, n) %>%
mutate(percent = round((yes / (yes + no)) * 100, 2))dots_std %>%
filter(session == "all") %>%
group_by(outlier, test, mask) %>%
summarize(n = n()) %>%
spread(outlier, n) %>%
mutate(percent = round((yes / (yes + no)) * 100, 2))# remove outlier column
betas_std = betas_std %>%
select(-outlier)
dots_std = dots_std %>%
select(-outlier)dataset = full_join(betas_std, dots_std, by = c("subjectID", "con", "process", "condition", "control", "session")) %>%
mutate(subjectID = as.character(subjectID)) %>%
left_join(., ind_diffs, by = "subjectID") %>%
ungroup() %>%
mutate(subjectID = as.factor(subjectID),
condition = as.factor(condition),
control = as.factor(control),
roi = as.factor(roi),
process = as.factor(process),
test = as.factor(test))
saveRDS(dataset, "full_dataset.RDS")
saveRDS(betas_std, "betas_dataset.RDS")
saveRDS(dots_std, "dots_dataset.RDS")dots_plot = dots_std %>%
filter(mask == "unmasked") %>%
filter(session == "all") %>%
select(-c(map, con, mask, session, dotProduct)) %>%
left_join(., ind_diffs) %>%
select(-c(sample, DBIC_ID, age, gender))
betas_plot = betas_std %>%
group_by(subjectID, process, condition, control) %>%
mutate(meanProcessPEstd = mean(meanPE_std, na.rm = TRUE)) %>%
select(-c(con, xyz, roi, meanPE, sdPE, meanPE_std)) %>%
unique() %>%
spread(process, meanProcessPEstd) %>%
mutate(balance = cognitive_control - reward) %>%
gather(process, meanProcessPEstd, cognitive_control, reward, value, balance) %>%
filter(session == "all") %>%
select(-c(session)) %>%
left_join(., ind_diffs) %>%
select(-c(sample, DBIC_ID, age, gender)) %>%
ungroup()
combined_plot = dots_plot %>%
rename("value_std" = dotProduct_std) %>%
mutate(test = ifelse(process == "craving_regulation", "neural signature", test)) %>%
select(subjectID, process, test, condition, control, value_std, fat, bmi) %>%
bind_rows(betas_plot %>%
rename("value_std" = meanProcessPEstd) %>%
mutate(test = "ROI") %>%
select(-contains("_z")))ema_tidy = ema %>%
select(subjectID, desire_pres0.prop, desire_str.m, goal_conf.m, resist.m, enact_prop, amount_Recoded.m, hunger.m) %>%
unique() %>%
mutate(enact_prop_z = as.numeric(scale(enact_prop)))
model_data = dataset %>%
filter(test == "association" & mask == "unmasked" & session == "all" & control == "nature" & condition == "food") %>%
unique() %>%
group_by(process, subjectID) %>%
mutate(meanProcessPEstd = mean(meanPE_std, na.rm = TRUE),
processPE = paste0(process, "_ROI"),
processPEV = paste0(process, "_PEV")) %>%
ungroup() %>%
select(-c(xyz, meanPE, meanPE_std, sdPE, dotProduct, map, process)) %>%
spread(processPEV, dotProduct_std) %>%
spread(processPE, meanProcessPEstd) %>%
group_by(subjectID) %>%
fill(contains("PEV"), .direction = "down") %>%
fill(contains("PEV"), .direction = "up") %>%
fill(contains("ROI"), .direction = "down") %>%
fill(contains("ROI"), .direction = "up") %>%
select(-c(roi, craving_ROI, craving_regulation_ROI, test, age, gender, mask)) %>%
unique() %>%
mutate(balance_ROI = cognitive_control_ROI - reward_ROI)
model_data_bmi = model_data %>%
select(-contains("fat")) %>%
na.omit()
model_data_fat = model_data %>%
select(-contains("bmi")) %>%
na.omit()
model_data_ema = model_data %>%
select(-contains("fat"), -contains("bmi")) %>%
right_join(., ema_tidy, by = "subjectID") %>%
na.omit()# source raincloud plot
source("https://gist.githubusercontent.com/benmarwick/2a1bb0133ff568cbe28d/raw/fb53bd97121f7f9ce947837ef1a4c65a73bffb3f/geom_flat_violin.R")
# plot
model_data %>%
left_join(., ema_tidy, by = "subjectID") %>%
filter(!is.na(subjectID)) %>%
select(subjectID, bmi, fat, enact_prop) %>%
rename("BMI" = bmi,
"body fat percentage" = fat,
"enactment" = enact_prop) %>%
gather(var, val, -subjectID) %>%
ggplot(aes("", val, fill = var)) +
geom_flat_violin(position = position_nudge(x = .1, y = 0), color = FALSE) +
geom_boxplot(width = .1, outlier.shape = NA) +
geom_point(position = position_jitter(width = .05), size = .5, alpha = .5) +
facet_wrap(~var, scales = "free") +
scale_fill_manual(values = pal_outcome) +
scale_x_discrete(expand = c(.1, 0)) +
coord_flip() +
labs(y = "\nvalue", x = "") +
theme_minimal(base_size = 14) +
theme(legend.position = "none")combined_plot %>%
filter(!process == "balance") %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature")),
process = gsub("_", " ", process),
test = factor(test, levels = c("neural signature", "ROI", "association", "uniformity"))) %>%
ggplot(aes(control, value_std, fill = condition)) +
stat_summary(fun.data = "mean_cl_boot", geom = "bar", position = position_dodge(.9)) +
stat_summary(fun.data = "mean_cl_boot", geom = "errorbar", position = position_dodge(.9), width = 0) +
facet_wrap(~ test + process, nrow = 3) +
scale_fill_manual(name = "", values = pal_condition) +
labs(x = "\ncontrol condition", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")combined_plot %>%
filter(!process == "balance") %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature")),
process = gsub("_", " ", process),
process = ifelse(!test == "ROI", sprintf("%s PEV", process), sprintf("%s ROI", process)),
test = ifelse(!test %in% c("association", "uniformity"), "none", test),
test = factor(test, levels = c("none", "association", "uniformity"))) %>%
ggplot(aes(control, value_std, fill = process)) +
stat_summary(fun.data = "mean_cl_boot", geom = "bar", position = position_dodge(.9)) +
stat_summary(fun.data = "mean_cl_boot", geom = "errorbar", position = position_dodge(.9), width = 0) +
facet_grid(test~condition) +
scale_fill_manual(name = "", values = pal_condition8) +
labs(x = "\ncontrol condition", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")Mean parameter estimates as a function of process, stimulus type, and control condition
betas_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(condition, meanProcessPEstd, fill = control)) +
stat_summary(fun.data = "mean_cl_boot", geom = "bar", position = position_dodge(.9)) +
stat_summary(fun.data = "mean_cl_boot", geom = "errorbar", position = position_dodge(.9), width = 0) +
facet_grid(~process) +
scale_fill_manual(name = "", values = pal_control) +
labs(x = "\nstimulus", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")betas_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(control, meanProcessPEstd, fill = condition)) +
stat_summary(fun.data = "mean_cl_boot", geom = "bar", position = position_dodge(.9)) +
stat_summary(fun.data = "mean_cl_boot", geom = "errorbar", position = position_dodge(.9), width = 0) +
facet_grid(~process) +
scale_fill_manual(name = "", values = pal_condition) +
labs(x = "\ncontrol condition", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")Correlations between variables and BMI or % fat
betas_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(bmi, meanProcessPEstd, color = condition)) +
geom_point(alpha = .2) +
geom_smooth(method = "lm", alpha = .25) +
facet_grid(control~process) +
scale_color_manual(name = "", values = pal_condition) +
labs(x = "\nBMI", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")betas_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(fat, meanProcessPEstd, color = condition)) +
geom_point(alpha = .2) +
geom_smooth(method = "lm", alpha = .25) +
facet_grid(control~process) +
scale_color_manual(name = "", values = pal_condition) +
labs(x = "\nbody composition (% fat)", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")Mean correspondance as a function of process, stimulus type, control condition, and test
dots_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(condition, dotProduct_std, fill = control)) +
stat_summary(fun.data = "mean_cl_boot", geom = "bar", position = position_dodge(.9)) +
stat_summary(fun.data = "mean_cl_boot", geom = "errorbar", position = position_dodge(.9), width = 0) +
facet_grid(test~process) +
scale_fill_manual(name = "", values = pal_control) +
labs(x = "\nstimulus", y = "mean correspondance\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")dots_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(control, dotProduct_std, fill = condition)) +
stat_summary(fun.data = "mean_cl_boot", geom = "bar", position = position_dodge(.9)) +
stat_summary(fun.data = "mean_cl_boot", geom = "errorbar", position = position_dodge(.9), width = 0) +
facet_grid(test~process) +
scale_fill_manual(name = "", values = pal_condition) +
labs(x = "\ncontrol condition", y = "mean correspondance\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")Correlations between variables and BMI or % fat
dots_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(bmi, dotProduct_std, color = condition)) +
geom_point(alpha = .2) +
geom_smooth(method = "lm", alpha = .25) +
facet_grid(control + test ~ process) +
scale_color_manual(name = "", values = pal_condition) +
labs(x = "\nBMI", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")dots_plot %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(fat, dotProduct_std, color = condition)) +
geom_point(alpha = .2) +
geom_smooth(method = "lm", alpha = .25) +
facet_grid(control + test ~ process) +
scale_color_manual(name = "", values = pal_condition) +
labs(x = "\nbody composition (% fat)", y = "mean standardized parameter estimate\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")model_data_ema %>%
select(subjectID, contains("ROI"), desire_pres0.prop, desire_str.m, goal_conf.m, resist.m, enact_prop, amount_Recoded.m) %>%
unique() %>%
gather(neuro_var, neuro_val, contains("ROI")) %>%
mutate(neuro_var = factor(neuro_var, levels = c("balance_ROI", "cognitive_control_ROI",
"reward_ROI", "value_ROI"))) %>%
gather(ema_var, ema_val, contains("m"), contains("prop")) %>%
group_by(ema_var) %>%
do({
plot = ggplot(., aes(neuro_val, ema_val, color = neuro_var)) +
geom_point(alpha = .5) +
geom_smooth(method = "lm", alpha = .1, fullrange = TRUE) +
scale_color_manual(name = "", values = pal_condition) +
labs(x = "\n mean parameter estimate", y = paste0(.$ema_var[[1]], "\n")) +
theme_minimal(base_size = 14)
print(plot)
data.frame()
})model_data_ema %>%
select(subjectID, contains("PEV"), desire_pres0.prop, desire_str.m, goal_conf.m, resist.m, enact_prop, amount_Recoded.m) %>%
unique() %>%
gather(neuro_var, neuro_val, contains("PEV")) %>%
mutate(neuro_var = factor(neuro_var, levels = c("cognitive_control_PEV", "reward_PEV", "value_PEV",
"craving_PEV", "craving_regulation_PEV"))) %>%
gather(ema_var, ema_val, contains("m"), contains("prop")) %>%
group_by(ema_var) %>%
do({
plot = ggplot(., aes(neuro_val, ema_val, color = neuro_var)) +
geom_point(alpha = .5) +
geom_smooth(method = "lm", alpha = .1, fullrange = TRUE) +
scale_color_manual(name = "", values = pal_condition[2:6]) +
labs(x = "\n mean parameter estimate", y = paste0(.$ema_var[[1]], "\n")) +
theme_minimal(base_size = 14)
print(plot)
data.frame()
})dots_plot %>%
filter(!grepl("craving", process)) %>%
left_join(., betas_plot) %>%
mutate(condition = factor(condition, levels = c("food", "dessert", "meal", "snack", "social", "nature"))) %>%
ggplot(aes(meanProcessPEstd, dotProduct_std, color = condition)) +
geom_point(alpha = .2) +
geom_smooth(method = "lm", alpha = .25) +
facet_grid(control + test ~ process) +
scale_color_manual(name = "", values = pal_condition) +
labs(x = "\nmean standardized parameter estimate", y = "mean standardized PEV\n") +
theme_minimal(base_size = 14) +
theme(axis.text.x = element_text(angle = 45, hjust = 1),
legend.position = "top",
legend.box = "horizontal")dots_plot %>%
unite(var, process, test, condition, control, sep = "_") %>%
mutate(var = sprintf("dot_%s", var)) %>%
spread(var, dotProduct_std) %>%
left_join(., betas_plot) %>%
unite(var, process, condition, control, sep = "_") %>%
mutate(var = sprintf("roi_%s", var)) %>%
spread(var, meanProcessPEstd) %>%
select(-subjectID) %>%
cor(., use = "pairwise.complete.obs") %>%
reshape2::melt() %>%
filter(!grepl("dot", Var1) & !grepl("roi", Var2)) %>%
ggplot(aes(Var1, Var2, fill = value)) +
geom_tile(color = "white") +
scale_fill_gradientn(name = "correlation\n", colors = c(pal_condition[1], "white", pal_condition[6]), limits = c(-1, 1), breaks = c(-1, 0, 1)) +
geom_text(aes(label = round(value, 1)), size = 2) +
labs(x = "", y = "") +
theme(legend.text = element_text(size = 8),
legend.position = "top",
legend.title = element_blank(),
axis.text = element_text(color = "black"),
axis.text.x = element_text(angle = 45, vjust = 1, hjust = 1),
axis.line = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
panel.background = element_blank())dots_plot %>%
unite(var, process, test, condition, control, sep = "_") %>%
mutate(var = sprintf("dot_%s", var)) %>%
spread(var, dotProduct_std) %>%
left_join(., betas_plot) %>%
unite(var, process, condition, control, sep = "_") %>%
mutate(var = sprintf("roi_%s", var)) %>%
spread(var, meanProcessPEstd) %>%
select(-subjectID) %>%
cor(., use = "pairwise.complete.obs") %>%
reshape2::melt() %>%
filter(grepl("roi", Var1) & grepl("roi", Var2)) %>%
ggplot(aes(Var1, Var2, fill = value)) +
geom_tile(color = "white") +
scale_fill_gradientn(name = "correlation\n", colors = c(pal_condition[1], "white", pal_condition[6]), limits = c(-1, 1), breaks = c(-1, 0, 1)) +
geom_text(aes(label = round(value, 1)), size = 2) +
labs(x = "", y = "") +
theme(legend.text = element_text(size = 8),
legend.position = "top",
legend.title = element_blank(),
axis.text = element_text(color = "black"),
axis.text.x = element_text(angle = 45, vjust = 1, hjust = 1),
axis.line = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
panel.background = element_blank())dots_plot %>%
unite(var, process, test, condition, control, sep = "_") %>%
mutate(var = sprintf("dot_%s", var)) %>%
spread(var, dotProduct_std) %>%
left_join(., betas_plot) %>%
unite(var, process, condition, control, sep = "_") %>%
mutate(var = sprintf("roi_%s", var)) %>%
spread(var, meanProcessPEstd) %>%
select(-subjectID) %>%
cor(., use = "pairwise.complete.obs") %>%
reshape2::melt() %>%
filter(grepl("dot", Var1) & grepl("dot", Var2)) %>%
ggplot(aes(Var1, Var2, fill = value)) +
geom_tile(color = "white") +
scale_fill_gradientn(name = "correlation\n", colors = c(pal_condition[1], "white", pal_condition[6]), limits = c(-1, 1), breaks = c(-1, 0, 1)) +
geom_text(aes(label = round(value, 1)), size = 2) +
labs(x = "", y = "") +
theme(legend.text = element_text(size = 8),
legend.position = "top",
legend.title = element_blank(),
axis.text = element_text(color = "black"),
axis.text.x = element_text(angle = 45, vjust = 1, hjust = 1),
axis.line = element_blank(),
panel.grid.major = element_blank(),
panel.grid.minor = element_blank(),
panel.border = element_blank(),
panel.background = element_blank())model_data %>%
left_join(., model_data_ema) %>%
gather(Variable, val, contains("PEV"), contains("ROI"), bmi, fat, enact_prop) %>%
filter(!Variable == "balance_ROI") %>%
group_by(Variable) %>%
summarize(M = mean(val, na.rm = TRUE),
SD = sd(val, na.rm = TRUE)) %>%
filter(!is.na(Variable)) %>%
mutate(Variable = gsub("_", " ", Variable),
Variable = ifelse(Variable == "bmi", "BMI",
ifelse(Variable == "fat", "Body fat percentage",
ifelse(Variable == "enact prop", "Enactment", Hmisc::capitalize(Variable)))),
Variable = factor(Variable, levels = c("BMI", "Body fat percentage", "Enactment",
"Reward ROI", "Reward PEV",
"Value ROI", "Value PEV",
"Craving regulation PEV", "Craving PEV",
"Cognitive control ROI", "Cognitive control PEV"))) %>%
arrange(Variable) %>%
kable(format = "pandoc", digits = 2)| Variable | M | SD |
|---|---|---|
| BMI | 23.19 | 3.36 |
| Body fat percentage | 25.43 | 8.00 |
| Enactment | 0.47 | 0.22 |
| Reward ROI | 0.06 | 0.50 |
| Reward PEV | 0.07 | 0.77 |
| Value ROI | 0.40 | 0.68 |
| Value PEV | 0.28 | 0.77 |
| Craving regulation PEV | 0.40 | 0.67 |
| Craving PEV | -0.13 | 0.76 |
| Cognitive control ROI | -0.13 | 0.50 |
| Cognitive control PEV | -0.11 | 0.78 |
cor_mat = model_data %>%
left_join(., model_data_ema) %>%
select(bmi, fat, enact_prop, contains("ROI"), contains("PEV"), -balance_ROI) %>%
gather(var, val, contains("PEV"), contains("ROI"), bmi, fat, enact_prop) %>%
mutate(var = gsub("_", " ", var),
var = ifelse(var == "bmi", "BMI",
ifelse(var == "fat", "Body fat percentage",
ifelse(var == "enact prop", "Enactment", Hmisc::capitalize(var)))),
var = factor(var, levels = c("BMI", "Body fat percentage", "Enactment",
"Reward ROI", "Reward PEV",
"Value ROI", "Value PEV",
"Craving regulation PEV", "Craving PEV",
"Cognitive control ROI", "Cognitive control PEV"))) %>%
arrange(var) %>%
spread(var, val) %>%
ungroup() %>%
select(-subjectID) %>%
as.matrix()
cors = Hmisc::rcorr(cor_mat)
cors$r %>%
as.data.frame() %>%
kable(format = "pandoc", digits = 2)| BMI | Body fat percentage | Enactment | Reward ROI | Reward PEV | Value ROI | Value PEV | Craving regulation PEV | Craving PEV | Cognitive control ROI | Cognitive control PEV | |
|---|---|---|---|---|---|---|---|---|---|---|---|
| BMI | 1.00 | 0.66 | -0.13 | -0.12 | -0.02 | 0.06 | 0.00 | 0.10 | -0.02 | -0.02 | -0.06 |
| Body fat percentage | 0.66 | 1.00 | -0.20 | -0.07 | -0.07 | 0.01 | -0.04 | 0.02 | -0.07 | -0.03 | -0.06 |
| Enactment | -0.13 | -0.20 | 1.00 | 0.21 | 0.31 | 0.23 | 0.29 | 0.03 | 0.21 | 0.01 | 0.05 |
| Reward ROI | -0.12 | -0.07 | 0.21 | 1.00 | 0.62 | 0.30 | 0.60 | 0.08 | 0.47 | 0.31 | 0.41 |
| Reward PEV | -0.02 | -0.07 | 0.31 | 0.62 | 1.00 | 0.58 | 0.91 | 0.19 | 0.81 | 0.47 | 0.59 |
| Value ROI | 0.06 | 0.01 | 0.23 | 0.30 | 0.58 | 1.00 | 0.65 | 0.28 | 0.51 | 0.06 | 0.16 |
| Value PEV | 0.00 | -0.04 | 0.29 | 0.60 | 0.91 | 0.65 | 1.00 | 0.15 | 0.76 | 0.30 | 0.40 |
| Craving regulation PEV | 0.10 | 0.02 | 0.03 | 0.08 | 0.19 | 0.28 | 0.15 | 1.00 | 0.24 | 0.21 | 0.23 |
| Craving PEV | -0.02 | -0.07 | 0.21 | 0.47 | 0.81 | 0.51 | 0.76 | 0.24 | 1.00 | 0.45 | 0.54 |
| Cognitive control ROI | -0.02 | -0.03 | 0.01 | 0.31 | 0.47 | 0.06 | 0.30 | 0.21 | 0.45 | 1.00 | 0.91 |
| Cognitive control PEV | -0.06 | -0.06 | 0.05 | 0.41 | 0.59 | 0.16 | 0.40 | 0.23 | 0.54 | 0.91 | 1.00 |
cors$P %>%
as.data.frame() %>%
mutate_if(is.numeric, funs(ifelse(between(., .05, .1), "†",
ifelse(between(., .01, .05), "*",
ifelse(between(., .001, .01), "**",
ifelse(. < .001, "***", "")))))) %>%
kable(format = "pandoc", digits = 2)| BMI | Body fat percentage | Enactment | Reward ROI | Reward PEV | Value ROI | Value PEV | Craving regulation PEV | Craving PEV | Cognitive control ROI | Cognitive control PEV |
|---|---|---|---|---|---|---|---|---|---|---|
| NA | *** | † | ||||||||
| *** | NA | † | ||||||||
| † | NA | * | ** | * | ** | * | ||||
| † | * | NA | *** | *** | *** | *** | *** | *** | ||
| ** | *** | NA | *** | *** | ** | *** | *** | *** | ||
| * | *** | *** | NA | *** | *** | *** | ** | |||
| ** | *** | *** | *** | NA | * | *** | *** | *** | ||
| ** | *** | * | NA | *** | *** | *** | ||||
| * | *** | *** | *** | *** | *** | NA | *** | *** | ||
| *** | *** | *** | *** | *** | NA | *** | ||||
| *** | *** | ** | *** | *** | *** | *** | NA |
model_data %>%
select(-contains("z"), -balance_ROI) %>%
left_join(., select(model_data_ema, subjectID, enact_prop)) %>%
ungroup() %>%
select_if(is.numeric) %>%
map_df(~ broom::tidy(t.test(.x, mu = 0, na.action = "na.omit")), .id = "x") %>%
rename("Variable" = x,
"t" = statistic,
"df" = parameter,
"p" = p.value) %>%
mutate(`Mdiff [95% CI]` = sprintf("%.02f [%.02f, %.02f]", estimate, conf.low, conf.high),
p = ifelse(p < .001, "< .001", gsub("0.(.*)", ".\\1", sprintf("%.3f", p))),
Variable = gsub("enact_prop", "Enactment", Variable),
Variable = gsub("fat", "Body fat percentage", Variable),
Variable = gsub("bmi", "BMI", Variable),
Variable = gsub("_", " ", Variable),
Variable = factor(Variable, levels = c("BMI", "Body fat percentage", "Enactment",
"reward ROI", "reward PEV",
"value ROI", "value PEV",
"craving regulation PEV", "craving PEV",
"cognitive control ROI", "cognitive control PEV"))) %>%
select(Variable, `Mdiff [95% CI]`, t, df, p) %>%
arrange(Variable) %>%
kable(format = "pandoc", digits = 2)| Variable | Mdiff [95% CI] | t | df | p |
|---|---|---|---|---|
| BMI | 23.19 [22.77, 23.60] | 109.65 | 252 | < .001 |
| Body fat percentage | 25.43 [24.44, 26.43] | 50.28 | 249 | < .001 |
| Enactment | 0.47 [0.42, 0.51] | 20.58 | 96 | < .001 |
| reward ROI | 0.06 [-0.00, 0.12] | 1.90 | 261 | .059 |
| reward PEV | 0.07 [-0.02, 0.17] | 1.51 | 261 | .132 |
| value ROI | 0.40 [0.31, 0.48] | 9.38 | 261 | < .001 |
| value PEV | 0.28 [0.18, 0.37] | 5.84 | 261 | < .001 |
| craving regulation PEV | 0.40 [0.31, 0.48] | 9.52 | 261 | < .001 |
| craving PEV | -0.13 [-0.22, -0.04] | -2.74 | 261 | .007 |
| cognitive control ROI | -0.13 [-0.19, -0.07] | -4.14 | 261 | < .001 |
| cognitive control PEV | -0.11 [-0.20, -0.01] | -2.20 | 261 | .029 |
options(na.action = "na.fail")
data_ctrl = caret::trainControl(method = "repeatedcv", number = 5, repeats = 5)bmi_roi_full = lm(bmi_z ~ value_ROI + reward_ROI + cognitive_control_ROI, data = model_data_bmi)
(bmi_roi_mods = MuMIn::dredge(bmi_roi_full))MuMIn::get.models(bmi_roi_mods, subset = TRUE) %>%
tibble() %>%
rename("model" = ".") %>%
mutate(tidied = purrr::map(model, broom::tidy),
model_num = row_number()) %>%
select(model_num, tidied) %>%
unnest()best_bmi_roi = caret::train(bmi_z ~ reward_ROI + value_ROI,
data = model_data_bmi,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)bmi_pev_full = lm(bmi_z ~ value_PEV + reward_PEV + cognitive_control_PEV + craving_PEV + craving_regulation_PEV, data = model_data_bmi)
(bmi_pev_mods = MuMIn::dredge(bmi_pev_full))MuMIn::get.models(bmi_pev_mods, subset = TRUE) %>%
tibble() %>%
rename("model" = ".") %>%
mutate(tidied = purrr::map(model, broom::tidy),
model_num = row_number()) %>%
select(model_num, tidied) %>%
unnest()best_bmi_pev = caret::train(bmi_z ~ craving_regulation_PEV,
data = model_data_bmi,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)bmi_null = lm(bmi_z ~ 1, data = model_data_bmi)
summary(bmi_null)##
## Call:
## lm(formula = bmi_z ~ 1, data = model_data_bmi)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.2100 -0.6542 -0.0558 0.4827 4.6115
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.01306 0.06326 0.206 0.837
##
## Residual standard error: 1.006 on 252 degrees of freedom
bmi_combined = caret::train(bmi_z ~ reward_ROI + value_ROI + craving_regulation_PEV,
data = model_data_bmi,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)
(bmi_aic = AIC(bmi_null, best_bmi_roi$finalModel, best_bmi_pev$finalModel, bmi_combined$finalModel) %>%
rownames_to_column() %>%
extract(rowname, "model", ".*bmi_(.*)\\$.*") %>%
mutate(model = toupper(model),
model = ifelse(is.na(model), "Null", model)) %>%
arrange(AIC))# ROI
summary(best_bmi_roi$finalModel)##
## Call:
## lm(formula = .outcome ~ ., data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1355 -0.6606 -0.0899 0.4662 4.5052
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.03113 0.07236 -0.430 0.6674
## reward_ROI -0.30523 0.12990 -2.350 0.0196 *
## value_ROI 0.15583 0.09654 1.614 0.1077
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9976 on 250 degrees of freedom
## Multiple R-squared: 0.02502, Adjusted R-squared: 0.01722
## F-statistic: 3.207 on 2 and 250 DF, p-value: 0.04214
best_bmi_roi## Linear Regression
##
## 253 samples
## 2 predictor
##
## No pre-processing
## Resampling: Cross-Validated (5 fold, repeated 5 times)
## Summary of sample sizes: 202, 201, 203, 203, 203, 203, ...
## Resampling results:
##
## RMSE Rsquared MAE
## 0.9949716 0.03211436 0.7469028
##
## Tuning parameter 'intercept' was held constant at a value of TRUE
sd(best_bmi_roi$resample$Rsquared)## [1] 0.03262441
sd(best_bmi_roi$resample$RMSE)## [1] 0.1099513
# combined
summary(bmi_combined$finalModel)##
## Call:
## lm(formula = .outcome ~ ., data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.1994 -0.6654 -0.0886 0.5098 4.4311
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.07322 0.07768 -0.942 0.3469
## reward_ROI -0.30708 0.12961 -2.369 0.0186 *
## value_ROI 0.11695 0.09989 1.171 0.2428
## craving_regulation_PEV 0.14370 0.09792 1.467 0.1435
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9953 on 249 degrees of freedom
## Multiple R-squared: 0.03338, Adjusted R-squared: 0.02173
## F-statistic: 2.866 on 3 and 249 DF, p-value: 0.03725
bmi_combined## Linear Regression
##
## 253 samples
## 3 predictor
##
## No pre-processing
## Resampling: Cross-Validated (5 fold, repeated 5 times)
## Summary of sample sizes: 200, 203, 204, 203, 202, 204, ...
## Resampling results:
##
## RMSE Rsquared MAE
## 0.9955154 0.03419623 0.7523059
##
## Tuning parameter 'intercept' was held constant at a value of TRUE
sd(bmi_combined$resample$Rsquared)## [1] 0.03657507
sd(bmi_combined$resample$RMSE)## [1] 0.1027421
fat_roi_full = lm(fat_z ~ value_ROI + reward_ROI + cognitive_control_ROI, data = model_data_fat)
(fat_roi_mods = MuMIn::dredge(fat_roi_full))MuMIn::get.models(fat_roi_mods, subset = TRUE) %>%
tibble() %>%
rename("model" = ".") %>%
mutate(tidied = purrr::map(model, broom::tidy),
model_num = row_number()) %>%
select(model_num, tidied) %>%
unnest()best_fat_roi = caret::train(fat_z ~ reward_ROI,
data = model_data_fat,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)fat_pev_full = lm(fat_z ~ value_PEV + reward_PEV + cognitive_control_PEV + craving_PEV + craving_regulation_PEV, data = model_data_fat)
(fat_pev_mods = MuMIn::dredge(fat_pev_full))MuMIn::get.models(fat_pev_mods, subset = TRUE) %>%
tibble() %>%
rename("model" = ".") %>%
mutate(tidied = purrr::map(model, broom::tidy),
model_num = row_number()) %>%
select(model_num, tidied) %>%
unnest()best_fat_pev = caret::train(fat_z ~ craving_PEV,
data = model_data_fat,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)fat_null = lm(fat_z ~ 1, data = model_data_fat)
fat_combined = caret::train(fat_z ~ reward_ROI + craving_PEV,
data = model_data_fat,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)
(fat_aic = AIC(fat_null, best_fat_roi$finalModel, best_fat_pev$finalModel, fat_combined$finalModel) %>%
rownames_to_column() %>%
extract(rowname, "model", ".*fat_(.*)\\$.*") %>%
mutate(model = toupper(model),
model = ifelse(is.na(model), "Null", model)) %>%
arrange(AIC))summary(fat_null)##
## Call:
## lm(formula = fat_z ~ 1, data = model_data_fat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.9531 -0.5945 0.1396 0.6457 2.4201
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 0.01679 0.06321 0.266 0.791
##
## Residual standard error: 0.9994 on 249 degrees of freedom
enact_prop_roi_full = lm(enact_prop_z ~ value_ROI + reward_ROI + cognitive_control_ROI, data = model_data_ema)
(enact_prop_roi_mods = MuMIn::dredge(enact_prop_roi_full))MuMIn::get.models(enact_prop_roi_mods, subset = TRUE) %>%
tibble() %>%
rename("model" = ".") %>%
mutate(tidied = purrr::map(model, broom::tidy),
model_num = row_number()) %>%
select(model_num, tidied) %>%
unnest()best_enact_prop_roi = caret::train(enact_prop_z ~ value_ROI,
data = model_data_ema,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)enact_prop_pev_full = lm(enact_prop_z ~ value_PEV + reward_PEV + cognitive_control_PEV + craving_PEV + craving_regulation_PEV, data = model_data_ema)
(enact_prop_pev_mods = MuMIn::dredge(enact_prop_pev_full))MuMIn::get.models(enact_prop_pev_mods, subset = TRUE) %>%
tibble() %>%
rename("model" = ".") %>%
mutate(tidied = purrr::map(model, broom::tidy),
model_num = row_number()) %>%
select(model_num, tidied) %>%
unnest()best_enact_prop_pev = caret::train(enact_prop_z ~ reward_PEV,
data = model_data_ema,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)enact_prop_null = lm(enact_prop_z ~ 1, data = model_data_ema)
enact_prop_combined = caret::train(enact_prop_z ~ value_ROI + reward_PEV,
data = model_data_ema,
trControl = data_ctrl,
method = "lm",
na.action = na.omit)
(enact_aic = AIC(enact_prop_null, best_enact_prop_roi$finalModel, best_enact_prop_pev$finalModel, enact_prop_combined$finalModel) %>%
rownames_to_column() %>%
extract(rowname, "model", ".*prop_(.*)\\$.*") %>%
mutate(model = toupper(model),
model = ifelse(is.na(model), "Null", model)) %>%
arrange(AIC))summary(best_enact_prop_pev$finalModel)##
## Call:
## lm(formula = .outcome ~ ., data = dat)
##
## Residuals:
## Min 1Q Median 3Q Max
## -2.17754 -0.58154 -0.01067 0.53885 2.40418
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -0.09627 0.09985 -0.964 0.33745
## reward_PEV 0.37158 0.11641 3.192 0.00192 **
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 0.9697 on 95 degrees of freedom
## Multiple R-squared: 0.09686, Adjusted R-squared: 0.08735
## F-statistic: 10.19 on 1 and 95 DF, p-value: 0.001917
best_enact_prop_pev## Linear Regression
##
## 97 samples
## 1 predictor
##
## No pre-processing
## Resampling: Cross-Validated (5 fold, repeated 5 times)
## Summary of sample sizes: 77, 79, 77, 78, 77, 78, ...
## Resampling results:
##
## RMSE Rsquared MAE
## 0.9722222 0.1498315 0.7677445
##
## Tuning parameter 'intercept' was held constant at a value of TRUE
sd(best_enact_prop_pev$resample$Rsquared)## [1] 0.1591469
sd(best_enact_prop_pev$resample$RMSE)## [1] 0.1075649
bmi_fit = data.frame(model = c("ROI", "PEV", "COMBINED"),
r2 = c(mean(best_bmi_roi$resample$Rsquared),
mean(best_bmi_pev$resample$Rsquared),
mean(bmi_combined$resample$Rsquared)),
r2_sd = c(sd(best_bmi_roi$resample$Rsquared),
sd(best_bmi_pev$resample$Rsquared),
sd(bmi_combined$resample$Rsquared)),
RMSE = c(mean(best_bmi_roi$resample$RMSE),
mean(best_bmi_pev$resample$RMSE),
mean(bmi_combined$resample$RMSE)),
RMSE_sd = c(sd(best_bmi_roi$resample$RMSE),
sd(best_bmi_pev$resample$RMSE),
sd(bmi_combined$resample$RMSE)))
bmi_table = broom::tidy(best_bmi_roi$finalModel, conf.int = TRUE) %>%
mutate(model = "ROI",
mod_sig = ifelse(pf(summary(best_bmi_roi$finalModel)$fstatistic[1],
summary(best_bmi_roi$finalModel)$fstatistic[2],
summary(best_bmi_roi$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", "")) %>%
bind_rows(broom::tidy(best_bmi_pev$finalModel, conf.int = TRUE) %>%
mutate(model = "PEV",
mod_sig = ifelse(pf(summary(best_bmi_pev$finalModel)$fstatistic[1],
summary(best_bmi_pev$finalModel)$fstatistic[2],
summary(best_bmi_pev$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", ""))) %>%
bind_rows(broom::tidy(bmi_combined$finalModel, conf.int = TRUE) %>%
mutate(model = "COMBINED",
mod_sig = ifelse(pf(summary(bmi_combined$finalModel)$fstatistic[1],
summary(bmi_combined$finalModel)$fstatistic[2],
summary(bmi_combined$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", ""))) %>%
bind_rows(broom::tidy(bmi_null, conf.int = TRUE) %>% mutate(model = "Null", mod_sig = "")) %>%
rename("SE" = std.error,
"t" = statistic,
"p" = p.value) %>%
mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", estimate, conf.low, conf.high),
p = ifelse(p < .001, "< .001", gsub("0.(.*)", ".\\1", sprintf("%.3f", p))),
term = gsub("\\(Intercept\\)", "Intercept", term),
term = gsub("_", " ", term),
term = Hmisc::capitalize(term)) %>%
left_join(., bmi_aic) %>%
left_join(., bmi_fit) %>%
mutate(`r2 (SD)` = ifelse(!model == "Null", sprintf("%.2f (%.2f)", r2, r2_sd), "--"),
`RMSE (SD)` = ifelse(!model == "Null", sprintf("%.2f (%.2f)", RMSE, RMSE_sd), "--")) %>%
unite(model, model, mod_sig, sep = "") %>%
select(model, AIC, `r2 (SD)`, `RMSE (SD)`, term, `b [95% CI]`, SE, t, p) %>%
arrange(AIC) %>%
mutate(model = gsub("COMBINED", "Combined", model),
outcome = "BMI") %>%
select(outcome, model, AIC, `r2 (SD)`, `RMSE (SD)`, term, `b [95% CI]`, SE, t, p)
bmi_table %>%
kable(format = "pandoc", digits = 2)| outcome | model | AIC | r2 (SD) | RMSE (SD) | term | b [95% CI] | SE | t | p |
|---|---|---|---|---|---|---|---|---|---|
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Intercept | -0.07 [-0.23, 0.08] | 0.08 | -0.94 | .347 |
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Reward ROI | -0.31 [-0.56, -0.05] | 0.13 | -2.37 | .019 |
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Value ROI | 0.12 [-0.08, 0.31] | 0.10 | 1.17 | .243 |
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Craving regulation PEV | 0.14 [-0.05, 0.34] | 0.10 | 1.47 | .144 |
| BMI | ROI* | 721.74 | 0.03 (0.03) | 0.99 (0.11) | Intercept | -0.03 [-0.17, 0.11] | 0.07 | -0.43 | .667 |
| BMI | ROI* | 721.74 | 0.03 (0.03) | 0.99 (0.11) | Reward ROI | -0.31 [-0.56, -0.05] | 0.13 | -2.35 | .020 |
| BMI | ROI* | 721.74 | 0.03 (0.03) | 0.99 (0.11) | Value ROI | 0.16 [-0.03, 0.35] | 0.10 | 1.61 | .108 |
| BMI | PEV | 723.45 | 0.03 (0.03) | 1.00 (0.13) | Intercept | -0.05 [-0.19, 0.10] | 0.07 | -0.67 | .507 |
| BMI | PEV | 723.45 | 0.03 (0.03) | 1.00 (0.13) | Craving regulation PEV | 0.16 [-0.03, 0.34] | 0.09 | 1.64 | .102 |
| BMI | Null | 724.15 | – | – | Intercept | 0.01 [-0.11, 0.14] | 0.06 | 0.21 | .837 |
fat_fit = data.frame(model = c("ROI", "PEV", "COMBINED"),
r2 = c(mean(best_fat_roi$resample$Rsquared),
mean(best_fat_pev$resample$Rsquared),
mean(fat_combined$resample$Rsquared)),
r2_sd = c(sd(best_fat_roi$resample$Rsquared),
sd(best_fat_pev$resample$Rsquared),
sd(fat_combined$resample$Rsquared)),
RMSE = c(mean(best_fat_roi$resample$RMSE),
mean(best_fat_pev$resample$RMSE),
mean(fat_combined$resample$RMSE)),
RMSE_sd = c(sd(best_fat_roi$resample$RMSE),
sd(best_fat_pev$resample$RMSE),
sd(fat_combined$resample$RMSE)))
fat_table = broom::tidy(best_fat_roi$finalModel, conf.int = TRUE) %>%
mutate(model = "ROI",
mod_sig = ifelse(pf(summary(best_fat_roi$finalModel)$fstatistic[1],
summary(best_fat_roi$finalModel)$fstatistic[2],
summary(best_fat_roi$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", "")) %>%
bind_rows(broom::tidy(best_fat_pev$finalModel, conf.int = TRUE) %>%
mutate(model = "PEV",
mod_sig = ifelse(pf(summary(best_fat_pev$finalModel)$fstatistic[1],
summary(best_fat_pev$finalModel)$fstatistic[2],
summary(best_fat_pev$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", ""))) %>%
bind_rows(broom::tidy(fat_combined$finalModel, conf.int = TRUE) %>%
mutate(model = "COMBINED",
mod_sig = ifelse(pf(summary(fat_combined$finalModel)$fstatistic[1],
summary(fat_combined$finalModel)$fstatistic[2],
summary(fat_combined$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", ""))) %>%
bind_rows(broom::tidy(fat_null, conf.int = TRUE) %>% mutate(model = "Null", mod_sig = "")) %>%
rename("SE" = std.error,
"t" = statistic,
"p" = p.value) %>%
mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", estimate, conf.low, conf.high),
p = ifelse(p < .001, "< .001", gsub("0.(.*)", ".\\1", sprintf("%.3f", p))),
term = gsub("\\(Intercept\\)", "Intercept", term),
term = gsub("_", " ", term),
term = Hmisc::capitalize(term)) %>%
left_join(., fat_aic) %>%
left_join(., fat_fit) %>%
mutate(`r2 (SD)` = ifelse(!model == "Null", sprintf("%.2f (%.2f)", r2, r2_sd), "--"),
`RMSE (SD)` = ifelse(!model == "Null", sprintf("%.2f (%.2f)", RMSE, RMSE_sd), "--")) %>%
unite(model, model, mod_sig, sep = "") %>%
select(model, AIC, `r2 (SD)`, `RMSE (SD)`, term, `b [95% CI]`, SE, t, p) %>%
arrange(AIC) %>%
mutate(model = gsub("COMBINED", "Combined", model),
outcome = "Body fat") %>%
select(outcome, model, AIC, `r2 (SD)`, `RMSE (SD)`, term, `b [95% CI]`, SE, t, p)
fat_table %>%
kable(format = "pandoc", digits = 2)| outcome | model | AIC | r2 (SD) | RMSE (SD) | term | b [95% CI] | SE | t | p |
|---|---|---|---|---|---|---|---|---|---|
| Body fat | Null | 712.18 | – | – | Intercept | 0.02 [-0.11, 0.14] | 0.06 | 0.27 | .791 |
| Body fat | ROI | 712.92 | 0.02 (0.02) | 1.00 (0.07) | Intercept | 0.02 [-0.10, 0.15] | 0.06 | 0.38 | .706 |
| Body fat | ROI | 712.92 | 0.02 (0.02) | 1.00 (0.07) | Reward ROI | -0.14 [-0.39, 0.11] | 0.13 | -1.12 | .265 |
| Body fat | PEV | 712.94 | 0.03 (0.03) | 1.00 (0.08) | Intercept | 0.00 [-0.12, 0.13] | 0.06 | 0.05 | .960 |
| Body fat | PEV | 712.94 | 0.03 (0.03) | 1.00 (0.08) | Craving PEV | -0.09 [-0.26, 0.07] | 0.08 | -1.11 | .268 |
| Body fat | Combined | 714.49 | 0.02 (0.03) | 1.01 (0.06) | Intercept | 0.01 [-0.12, 0.14] | 0.07 | 0.19 | .849 |
| Body fat | Combined | 714.49 | 0.02 (0.03) | 1.01 (0.06) | Reward ROI | -0.09 [-0.38, 0.19] | 0.14 | -0.66 | .508 |
| Body fat | Combined | 714.49 | 0.02 (0.03) | 1.01 (0.06) | Craving PEV | -0.06 [-0.25, 0.13] | 0.10 | -0.65 | .516 |
enact_prop_fit = data.frame(model = c("ROI", "PEV", "COMBINED"),
r2 = c(mean(best_enact_prop_roi$resample$Rsquared),
mean(best_enact_prop_pev$resample$Rsquared),
mean(enact_prop_combined$resample$Rsquared)),
r2_sd = c(sd(best_enact_prop_roi$resample$Rsquared),
sd(best_enact_prop_pev$resample$Rsquared),
sd(enact_prop_combined$resample$Rsquared)),
RMSE = c(mean(best_enact_prop_roi$resample$RMSE),
mean(best_enact_prop_pev$resample$RMSE),
mean(enact_prop_combined$resample$RMSE)),
RMSE_sd = c(sd(best_enact_prop_roi$resample$RMSE),
sd(best_enact_prop_pev$resample$RMSE),
sd(enact_prop_combined$resample$RMSE)))
enact_table = broom::tidy(best_enact_prop_roi$finalModel, conf.int = TRUE) %>%
mutate(model = "ROI",
mod_sig = ifelse(pf(summary(best_enact_prop_roi$finalModel)$fstatistic[1],
summary(best_enact_prop_roi$finalModel)$fstatistic[2],
summary(best_enact_prop_roi$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", "")) %>%
bind_rows(broom::tidy(best_enact_prop_pev$finalModel, conf.int = TRUE) %>%
mutate(model = "PEV",
mod_sig = ifelse(pf(summary(best_enact_prop_pev$finalModel)$fstatistic[1],
summary(best_enact_prop_pev$finalModel)$fstatistic[2],
summary(best_enact_prop_pev$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", ""))) %>%
bind_rows(broom::tidy(enact_prop_combined$finalModel, conf.int = TRUE) %>%
mutate(model = "COMBINED",
mod_sig = ifelse(pf(summary(enact_prop_combined$finalModel)$fstatistic[1],
summary(enact_prop_combined$finalModel)$fstatistic[2],
summary(enact_prop_combined$finalModel)$fstatistic[3], lower = FALSE) < .05, "*", ""))) %>%
bind_rows(broom::tidy(enact_prop_null, conf.int = TRUE) %>% mutate(model = "Null", mod_sig = "")) %>%
rename("SE" = std.error,
"t" = statistic,
"p" = p.value) %>%
mutate(`b [95% CI]` = sprintf("%.2f [%.2f, %.2f]", estimate, conf.low, conf.high),
p = ifelse(p < .001, "< .001", gsub("0.(.*)", ".\\1", sprintf("%.3f", p))),
term = gsub("\\(Intercept\\)", "Intercept", term),
term = gsub("_", " ", term),
term = Hmisc::capitalize(term)) %>%
left_join(., enact_aic) %>%
left_join(., enact_prop_fit) %>%
mutate(`r2 (SD)` = ifelse(!model == "Null", sprintf("%.2f (%.2f)", r2, r2_sd), "--"),
`RMSE (SD)` = ifelse(!model == "Null", sprintf("%.2f (%.2f)", RMSE, RMSE_sd), "--")) %>%
unite(model, model, mod_sig, sep = "") %>%
select(model, AIC, `r2 (SD)`, `RMSE (SD)`, term, `b [95% CI]`, SE, t, p) %>%
arrange(AIC) %>%
mutate(model = gsub("COMBINED", "Combined", model),
outcome = "Enactment") %>%
select(outcome, model, AIC, `r2 (SD)`, `RMSE (SD)`, term, `b [95% CI]`, SE, t, p)
enact_table %>%
kable(format = "pandoc", digits = 2)| outcome | model | AIC | r2 (SD) | RMSE (SD) | term | b [95% CI] | SE | t | p |
|---|---|---|---|---|---|---|---|---|---|
| Enactment | PEV* | 273.28 | 0.15 (0.16) | 0.97 (0.11) | Intercept | -0.10 [-0.29, 0.10] | 0.10 | -0.96 | .337 |
| Enactment | PEV* | 273.28 | 0.15 (0.16) | 0.97 (0.11) | Reward PEV | 0.37 [0.14, 0.60] | 0.12 | 3.19 | .002 |
| Enactment | Combined* | 274.60 | 0.12 (0.13) | 0.97 (0.13) | Intercept | -0.15 [-0.38, 0.09] | 0.12 | -1.25 | .215 |
| Enactment | Combined* | 274.60 | 0.12 (0.13) | 0.97 (0.13) | Value ROI | 0.13 [-0.19, 0.46] | 0.16 | 0.81 | .419 |
| Enactment | Combined* | 274.60 | 0.12 (0.13) | 0.97 (0.13) | Reward PEV | 0.32 [0.05, 0.58] | 0.14 | 2.33 | .022 |
| Enactment | ROI* | 278.06 | 0.10 (0.09) | 0.99 (0.14) | Intercept | -0.19 [-0.43, 0.05] | 0.12 | -1.59 | .115 |
| Enactment | ROI* | 278.06 | 0.10 (0.09) | 0.99 (0.14) | Value ROI | 0.33 [0.04, 0.61] | 0.14 | 2.27 | .026 |
| Enactment | Null | 281.16 | – | – | Intercept | -0.04 [-0.25, 0.16] | 0.10 | -0.42 | .676 |
bind_rows(bmi_table, fat_table, enact_table) %>%
kable(format = "pandoc", digits = 2)| outcome | model | AIC | r2 (SD) | RMSE (SD) | term | b [95% CI] | SE | t | p |
|---|---|---|---|---|---|---|---|---|---|
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Intercept | -0.07 [-0.23, 0.08] | 0.08 | -0.94 | .347 |
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Reward ROI | -0.31 [-0.56, -0.05] | 0.13 | -2.37 | .019 |
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Value ROI | 0.12 [-0.08, 0.31] | 0.10 | 1.17 | .243 |
| BMI | Combined* | 721.56 | 0.03 (0.04) | 1.00 (0.10) | Craving regulation PEV | 0.14 [-0.05, 0.34] | 0.10 | 1.47 | .144 |
| BMI | ROI* | 721.74 | 0.03 (0.03) | 0.99 (0.11) | Intercept | -0.03 [-0.17, 0.11] | 0.07 | -0.43 | .667 |
| BMI | ROI* | 721.74 | 0.03 (0.03) | 0.99 (0.11) | Reward ROI | -0.31 [-0.56, -0.05] | 0.13 | -2.35 | .020 |
| BMI | ROI* | 721.74 | 0.03 (0.03) | 0.99 (0.11) | Value ROI | 0.16 [-0.03, 0.35] | 0.10 | 1.61 | .108 |
| BMI | PEV | 723.45 | 0.03 (0.03) | 1.00 (0.13) | Intercept | -0.05 [-0.19, 0.10] | 0.07 | -0.67 | .507 |
| BMI | PEV | 723.45 | 0.03 (0.03) | 1.00 (0.13) | Craving regulation PEV | 0.16 [-0.03, 0.34] | 0.09 | 1.64 | .102 |
| BMI | Null | 724.15 | – | – | Intercept | 0.01 [-0.11, 0.14] | 0.06 | 0.21 | .837 |
| Body fat | Null | 712.18 | – | – | Intercept | 0.02 [-0.11, 0.14] | 0.06 | 0.27 | .791 |
| Body fat | ROI | 712.92 | 0.02 (0.02) | 1.00 (0.07) | Intercept | 0.02 [-0.10, 0.15] | 0.06 | 0.38 | .706 |
| Body fat | ROI | 712.92 | 0.02 (0.02) | 1.00 (0.07) | Reward ROI | -0.14 [-0.39, 0.11] | 0.13 | -1.12 | .265 |
| Body fat | PEV | 712.94 | 0.03 (0.03) | 1.00 (0.08) | Intercept | 0.00 [-0.12, 0.13] | 0.06 | 0.05 | .960 |
| Body fat | PEV | 712.94 | 0.03 (0.03) | 1.00 (0.08) | Craving PEV | -0.09 [-0.26, 0.07] | 0.08 | -1.11 | .268 |
| Body fat | Combined | 714.49 | 0.02 (0.03) | 1.01 (0.06) | Intercept | 0.01 [-0.12, 0.14] | 0.07 | 0.19 | .849 |
| Body fat | Combined | 714.49 | 0.02 (0.03) | 1.01 (0.06) | Reward ROI | -0.09 [-0.38, 0.19] | 0.14 | -0.66 | .508 |
| Body fat | Combined | 714.49 | 0.02 (0.03) | 1.01 (0.06) | Craving PEV | -0.06 [-0.25, 0.13] | 0.10 | -0.65 | .516 |
| Enactment | PEV* | 273.28 | 0.15 (0.16) | 0.97 (0.11) | Intercept | -0.10 [-0.29, 0.10] | 0.10 | -0.96 | .337 |
| Enactment | PEV* | 273.28 | 0.15 (0.16) | 0.97 (0.11) | Reward PEV | 0.37 [0.14, 0.60] | 0.12 | 3.19 | .002 |
| Enactment | Combined* | 274.60 | 0.12 (0.13) | 0.97 (0.13) | Intercept | -0.15 [-0.38, 0.09] | 0.12 | -1.25 | .215 |
| Enactment | Combined* | 274.60 | 0.12 (0.13) | 0.97 (0.13) | Value ROI | 0.13 [-0.19, 0.46] | 0.16 | 0.81 | .419 |
| Enactment | Combined* | 274.60 | 0.12 (0.13) | 0.97 (0.13) | Reward PEV | 0.32 [0.05, 0.58] | 0.14 | 2.33 | .022 |
| Enactment | ROI* | 278.06 | 0.10 (0.09) | 0.99 (0.14) | Intercept | -0.19 [-0.43, 0.05] | 0.12 | -1.59 | .115 |
| Enactment | ROI* | 278.06 | 0.10 (0.09) | 0.99 (0.14) | Value ROI | 0.33 [0.04, 0.61] | 0.14 | 2.27 | .026 |
| Enactment | Null | 281.16 | – | – | Intercept | -0.04 [-0.25, 0.16] | 0.10 | -0.42 | .676 |
# BMI
car::vif(bmi_combined$finalModel)## reward_ROI value_ROI craving_regulation_PEV
## 1.101738 1.185011 1.085275
# %fat
car::vif(fat_combined$finalModel)## reward_ROI craving_PEV
## 1.303633 1.303633
# enactment
car::vif(enact_prop_combined$finalModel)## value_ROI reward_PEV
## 1.348297 1.348297